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ABSTRACT 


In  this  paper  an  analytical  technique  based  on  the  method  of  undetermined 
coefficients  is  applied  to  the  problem  of  computing  the  theoretical  maximum 
entropy  (MEM)  spectral  estimate  when  the  correlation  function  of  the  data  is 
known  exactly  and  corresponds  to  N sinusoids  in  additive  white  noise  and  to  N 
sinusoids  in  additive  1-pole,  low  pass  noise.  For  the  white  noise  case,  the  L 
prediction  filter  coefficients  are  expanded  directly  in  terms  of  the  input  sinusoids. 
This  expansion  leads  to  a transformation  of  the  L X L normal  equations  for  the 
prediction  filter  coefficients  to  a set  of  2N  X 2N  equations.  The  transformed 
equations  are  a smaller  set  of  equations  to  be  solved  whenever  L > 2N  and  pro- 
vide a convenient  description  of  the  interaction  between  the  various  frequency 
components  of  the  sinusoids  which  occurs  in  the  MEM  estimate.  Fvirther,  for 
certain  cases  where  there  is  little  interaction  between  some  of  the  frequency 
components  of  the  sinusoids,  the  solution  of  the  2N  X 2N  equations  may  be 
approximated  (to  zeroth  order)  by  the  solution  of  a smaller  set  of  coupled  equa- 
tions. A better  approximation  to  the  exact  solution  of  the  2N  X 2N  equations 
can  then  be  obtained  from  a perturbation  expansion  of  the  exact  solution  about 
the  zeroth  order  approximation. 
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For  the  case  of  N sinusoids  in  1-pole,  low  pass  noise,  the  L prediction 
filter  coefficients  are  expanded  in  terms  of  the  input  sinusoids  as  well  as  two 


delta  functions  which  occur  at  the  beginning  and  end  of  the  filter.  This  expansion 
also  leads  to  a set  of  2N  X 2N  equations.  For  this  case  the  values  of  the  MEM 
estimate  evaluated  at  the  frequencies  of  the  sinusoids  are  shown  to  be  a function 
of  the  frequencies  of  the  sinusoids.  This  result  is  reasonable  since  the  signal- 
to-noise  ratio  per  unit  bandwidth  is  also  a function  of  frequency. 
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INTRODUCTION 


! 

t 

i 

r 

t 
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The  maximum  entropy  method  (MEM)  of  spectral  analysis,  which  was  origi- 
nally proposed  by  Burg  (1967,  1975),  has  been  widely  applied  in  geophysical  data 
processing.  The  correspondence  between  MEM  and  linear  prediction  filtering, 
as  discussed  by  van  den  Bos  (1971),  has  allowed  the  application  of  a large  body 
of  literature  on  autoregressive  (AR)  time  series  analysis  to  MEM.  Ulrych  and 
Bishop  (1975)  and  Ulrych  and  Clayton  (1976)  give  a thorough  review  and  discussion 
of  MEM  and  AR  analysis. 

An  important  Explication  of  MEM  is  to  the  spectrum  analysis  of  data  contain- 
ing multiple  sinusoids  in  noise.  Even  though  the  proper  time  series  model  for 
this  case  is  not  an  AR  model  (as  will  be  discussed  in  more  detail  in  the  next  sec- 
tions), MEM  can  still  provide  excellent  sinusoid  resolution  (especially  for  large 
signal-to-noise  ratios).  The  improved  resolution  offered  by  MEM  over  the  more 
conventional  Fourier  spectral  estimation  methods  for  this  cEise  has  been  well  docu- 
mented in  the  literature  (e.  g. , Ulrych,  1972;  Ulrych  and  Bishop,  1975).  For  pur- 
poses of  comparing  the  performance  of  MEM  with  other  spectral  estimation  tech- 
niques, Lacoss  (1971)  and  more  recently  Frost  (1977)  and  Marple  (1976)  have 
examined  the  theoretical  MEM  spectral  estimate  when  the  correlation  function 
of  the  data  is  known  exactly  and  consists  of  sinusoids  in  white  noise.  The  use  of 
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an  exact  matrix  inverse  identity  for  computing  the  theoretical  MEM  spectral 
estimate  for  the  cases  of  one  and  two  sinusoids  in  white  noise  was  discussed  by 
Lacoss  (1971)  and  Frost  (1977).  However,  when  there  are  more  than  two  sinus- 
oids, the  use  of  the  inverse  identity  becomes  tedious.  Even  for  two  sinusoids, 
the  analytical  form  of  the  theoretical  MEM  spectral  estimate  which  is  obtained 
from  the  repeated  application  of  the  inverse  identity  provides  little  insight  into 
the  interaction  between  the  positive  and  negative  frequency  components  of  the 
sinusoids  which  occurs  in  the  MEM  estimate  (Marple,  1976).  In  this  paper  an 
alternative  approach  based  on  the  method  of  undetermined  coefficients  will  be 
used  to  compute  the  theoretical  MEM  spectral  estimate  when  the  correlation 
lags  of  the  data  are  known  exactly  and  correspond  to  multiple  sinusoids  in  white 
noise  and  to  multiple  sinusoids  in  1-pole,  low  pass  noise.  It  will  be  shown  that 
the  method  of  undetermined  coefficients  provides  a convenient  description  of  the 
interaction  between  the  various  frequency  components  of  the  sinusoids  which  oc- 
curs in  the  theoretical  MEM  estimate.  Further,  useful  approximations  for  the 
theoretical  MEM  spectral  estimate  may  be  obtained  for  cases  in  which  the  inter- 
action between  certain  sinusoidal  frequency  components  is  small  but  not  negli- 
gible. For  the  case  of  1-pole,  low  pass  noise,  it  will  be  shown  that  the  MEM 
spectral  resolution  is  a function  of  the  frequency  and  becomes  worse  in  the 
spectral  regions  where  the  noise  power  increases. 


MEM  SPECTRAL  ANALYSIS  OF  MULTIPLE  SINUSOIDS  IN  WHITE  NOISE 


In  this  section  we  will  consider  the  case  when  the  correlation  lags  of  the 
data,  P(-^),  -2=0,1, . . . ,L,  are  known  exactly  and  can  be  written  as  follows 


I'! 

P(-2)  = (>{&)  + £ cos  2fff 

o n n 


where  6(-2)  is  the  Kronecker  delta  function,  is  the  power  in  the  n*^^  sinusoid, 

n ’ 

f^  represents  the  frequencies  of  the  sinusoids  (which  are  normalized  to  the  sample 
2 

frequency),  and  is  the  white  noise  power.  The  MEM  spectral  estimate, 

S(f),  may  be  written  as  follows  (tflrych  and  Bishop,  1975) 


S(f)  = Aa^  Q(f) 


where  A is  the  sampling  interval,  and 


k=0 


where,  f is  normalized  to  the  sample  frequency.  The  L prediction  coeffi- 
cients, g(k),k=0, . . . ,L-1,  are  obtained  from  the  normal  equations 


t 
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(4) 


L-1 

E p(X-k)  g(k)  = P(^+l),  1=0, . . . ,L-1 
k=0 


and  the  constant  O 

V 


is  obtained  from  the  g(k)  by 


2 

a = p(0)  - E g(k)  p(k+i)  (5) 

k=0 


Lacoss  (1971)  and  Marple  (1976)  treat  the  problem  of  computing  S(f)  either 
through  the  direct  numerical  solution  of  equations  (4)-(5)  or  through  the  use  of 
a well  known  matrix  inversion  identity  which  is  sometimes  referred  to  as  Wood- 
bury's identity  (Zielke,  1968).  As  noted  above,  the  application  of  Woodbury's 
identity  becomes  quite  tedious  and  leads  to  very  extensive  anal3d;ic  expressions 
for  g(k)  and  S(f)  if  N is  larger  than  2. 

An  alternative  technique  of  examining  S(f)  analytically  and  numerically  is 
the  method  of  undetermined  coefficients.  This  method  consists  of  substituting 
a solution  for  g(k),  which  is  expressed  in  terms  of  unknown  constants,  into 
equation  (4).  This  substitution  then  leads  to  a set  of  equations  for  the  unknown 
constants.  This  technique  was  originally  applied  to  the  continuous  analog  of 
equation  (4)  by  Zadeh  and  Ragazzini  (1950)  and  more  recently  has  been  applied 
directly  to  equation  (4)  by  Satorius  and  Zeidler  (1977)  when  the  spectral  density 
of  the  data  contains  both  poles  and  zeroes. 
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The  form  of  the  assumed  solution  for  g(k)  when  p(-2)  is  given  by  equation 


(1)  is 


271  jf  k 
n 


g(k)  = E A e 
n=l 


where  we  define  (for  notat’onal  convenience)  f ^ = -f  (n=l,  2, . . . ,N) 

n+N  n 

and  the  are  to  be  determined.  This  particular  choice  for  the  assumed 
solution  leads  to  precisely  2N  equations  for  the  A^.  Substituting  equation  (6) 
into  (4)  with  p(-^)  given  by  equation  (1)  and  equating  coefficients  of  exp(277jf^£) 

(for  r=l,2, . . . , 2N)  in  the  resulting  equation  leads  to  the  following  2N  equations 
for  the  A 


n=l  L+2a  /a 

, or 

nT^r 


2 2 

where  in  (7)  we  have  defined  ^ (n=l,2, . . . ,N)  and  y is  given  by 

' ' n+N  n m 


y = L n r 

" h*2aW 

o r 
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where 


L-1 

0 (f)  = E e 
^ k=0 


27rjfk 


(9) 


It  is  noted  that  the  net  effect  of  the  method  of  undetermined  coefficients  is  to 
transform  the  original  LXL  equations  (equation  (4))  to  the  set  of  2N  X 2N  equa- 
tions (equation  (7)).  This  transformation  yields  a smaller  set  of  equations  to  be 
solved  whenever  L>  2N. 

Equations  (6)-(8)  show  that  the  prediction  coefficients,  g(k),  can  be  ex- 
pressed as  a sum  of  the  positive  and  negative  frequency  components  of  the  input 
sinusoids  and  that  the  amplitude  of  each  sinusoid,  A^,  is  coupled  to  the  amplitude 
of  all  the  other  sinusoids  through  coupling  coefficients,  . The  coupling 

coefficients  vanish  if  ~ is  some  integral  multiple  of  1/L.  Also,  as 
2 2 

the  factor,  L + 20  /a  , becomes  large,  the  y approach  zero.  As  the  y 
u r rn  rn 

^proach  zero,  equation  (7)  decouples  and  the  A^  are  given  to  a good  approxi- 
mation by 


A 


r 


e 

2 , 2 
L+2a  /a 
o r 


r=l,2. 


2N. 


(10) 
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Using  equations  (6)-(8)  we  may  express  all  quantities  of  interest  in  terms 

2 

of  the  A . From  (5)  we  have  for  a 
n V 


2 2 2N  -217jf 

a^=a  i+SAe 
VO  , n 

n=l 


As  the  y approach  zero,  CT  is  given  to  a good  approximation  by 
m V 


2 2 / 2 , 2 


From  equation  (3),  we  have  for  Q(f) 


1 -2ffjf  . 1-2 

= ll-e  S A„»  ,t  -t)l  . 

n=l 


It  is  interesting  to  note  that  Q(f)  evaluated  at  the  frequency  of  the  r sinusoid 

2 2 

can  be  simply  expressed  in  terms  of  A and  O’  /o  . The  result  for  Q(f  ) is 

r o r r 


Q(V  = r=l N. 
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Equation  (14)  is  valid  regardless  of  the  frequency  separation  between  the  N 


sinusoids.  (Of  course,  as  the  frequency  separations  approach  zero,  the  maxima 
of  Q(f)  will  not  necessarily  be  equal  to  the  Q(fp.)  As  the  approach  zero, 
Q(fp  is  given  approximately  by 

2 

Q(y"‘  (t^V4aJ)(L  + 2aW)  , r=l N.  (15) 

2 2 

When  L » 2CT^/a^,  equation  (15)  becomes 

Q(y  r=l,...,N.  (16) 

A result  similar  to  equation  (16)  was  also  obtained  by  Lacoss  (1971)  for  the 
theoretical  peak  values  of  the  MEM  estimate  of  a real  sinusoid  in  white  noise. 

As  pointed  out  in  the  recent  paper  by  Ulrych  and  Clayton  (1976),  the  proper 

time  series  model  for  N sinusoids  in  white  noise  is  an  autoregressive-moving 

average  (ARMA)  model  which  contains  2N  autoregressive  (AR)  terms  and  2N 

moving  average  (MA)  terms.  Therefore,  since  an  infinite  order  AR  model  is 

required  to  model  a finite  order  ARMA  process  (Gersch  and  Sharpe,  1973), 

2 ^ . 

it  is  expected  that  as  L ®,  S(f)  will  converge  to  (for  f f f^;  n=] N). 

This  can  be  seen  quite  simply  by  noting  from  equation  (10)  that 


12 


provided  f (n-l,...,N).  Therefore,  from  equations  (2),  (11),  and  (13), 

S(f)  converges  pointwise  to  ^.s  L“*  “ when  f ^ f and  from  equation  (16) 

it  is  seen  that  provided  the  correlation  lags  are  known  exactly,  the  resolution 

capabilities  of  S(f)  inprove  without  limit  as  L-  “ regardless  of  the  value  of 
2 

^o*  course,  in  reality  the  correlation  lags  are  never  known  exactly  but  must 
be  estimated  from  a finite  amount  of  data  and,  therefore,  a practical  limit  is 
imposed  on  the  value  of  L which  provides  a tradeoff  between  MEM  resolution 
and  the  confidence  in  the  estimates  of  the  correlation  lags  (as  well  as  the  in- 
creased inaccuracies  which  occur  when  computing  a larger  number  of  predic- 
tion filter  coefficients).  The  criteria  which  have  been  most  frequently  applied 
to  the  problem  of  determining  the  optimum  value  of  L are  the  Akaike  final 
prediction  error  (FPE)  criterion  (Akaike,  1969,  1970)  and  the  Akaike  informa- 
tion theoretic  (AIC)  criterion  (Akaike,  1972,  1974).  The  application  of  these 
different  criteria  has  been  considered  by  a number  of  authors  (e.g. , Gersch  and 
Sharpe,  1973;  Akaike,  1974;  Jones,  1976).  The  problem  of  determining  the 
optimal  value  of  L is  further  discussed  by  Ulrych  and  Bishop  (1975)  and 
Ulrych  and  Clayton  (1976). 


13 


APPROXIMATIONS  TO  THE  THEORETICAL  MEM  SPECTRAL 
ESTIMATE  FOR  SINUSOIDS  IN  WHITE  NOISE 


Note  that  when  some  of  the  are  negligible,  an  approximate  solution 

(zeroth  order)  for  the  A may  be  obtained  by  setting  the  negligible  y to  zero 

m 

in  equation  (7).  A better  approximation  to  the  true  A^  may  be  obtained  by  ex- 
panding equation  (7)  in  a perturbation  expansion  about  the  zeroth  order  approxi- 
mation. In  particular,  let  B be  the  matrix  with  elements  B given  by 

kj 


1;  if  k = j 

\j’  "^^Sligible  in  zeroth  order 

0;  otherwise 


Further,  lety^^^ 


be  the  matrix  with  elements 


given  by 


0;  if  k = j 

if  is  negligible  in 
0;  otherwise 


zeroth  order 


Equation  (7)  may  now  be  expressed  in  the  equivalent  matrix  form 


(17) 
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I 


where  I is  the  2N  x 2N  identity  matrix;  A is  a column  vector  with 
element  g^ven  by  Aj^;  and  _A^  ^ = B . F,  where  F is  a column  vector  with 
k**’  element  given  by  exp(2irj^)/(L  + 2(7^/a^) . a<®>  is  the  zeroth  order 

approximation  to  A.  Equation  (17)  can  be  further  expanded  in  a perturbation 
expansion  for  A in  terms  of  the  matrix  The  result  is 


(18) 


(19) 


Equation  (18)  is  the  desired  perturbation  expansion  for  A in  powers  of 
The  convergence  of  equation  (18)  is  discussed  in  the  Appendix. 

As  a specific  example,  consider  the  case  of  two  sinusoids  in  white  noise  in 
which  there  is  little  interaction  between  the  positive  and  negative  frequency  com- 
ponents of  the  sinusoids  but  there  is  appreciable  interaction  between  the  two  posi- 
tive frequencies  (and,  therefore,  the  two  negative  frequencies).  For  this  case, 
the  8 coefficients  y^^,  y^^,  y^^,  y^^,  and  y^^  may  be  neglected 

in  zeroth  order,  and,  therefore,  A^  ^ is  obtained  by  solving  the  two  independent 

sets  of  2 X 2 equations  involving  only  the  non -negligible  coefficients  y , y , 

X2  21 
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y , and  y . From  equations  (18)-(19),  it  is  strai^tforward  to  obtain  approxi- 
mations  for  Q(f).  To  zeroth  order,  Q(f)  is  given  by 

Q(f)  ^ £ 1 1 - -f)  + 1 . (20) 

In  equation  (20),  we  have  neglected  A^^^  and  A^®^  which  is  equivalent  to  neglect- 
ing  the  negative  part  of  the  frequency  spectrum  in  Q(f)  and  is  consistent  with  the 
neglection  of  the  interaction  between  the  positive  and  negative  frequency  compo- 
nents in  the  zeroth  order  approximation.  The  approximations  in  obtaining  equa- 
tion (20)  are  equivalent  to  applying  Woodbury's  identity  to  the  correlation  matrix 
formed  by  the  positive  frequency  components  of  the  two  sinusoids  as  suggested 
by  Lacoss  (1971).  Equation  (20)  will  give  a good  approximation  to  Q(f)  for  f 
sufficiently  far  from  zero  and  will  also  give  the  zeroth  order  approximation  to 
Q(f^),  i.e. , 


(21) 


However,  for  certain  applications  (e.g. , Frost,  1977),  the  errors  introduced  by 
neglecting  the  negative  frequency  components  in  Q(f)  can  be  appreciable.  To 
include  these  effects  (to  first  order)  we  have  from  equations  (18)-(19) 


Q(f) 


= |l 


n=l 


-'•I 

\ n n / L n 


-2 


(22) 
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Equation  (22)  will  provide  a better  approximation  to  Q(f)  over  a wider  range  of 
frequencies  than  \f);  however,  will  only  be  correct  to  zeroth  order, 

i.e.,  (from  (19)  and  (22)) 


Note  that  for  small  values  of  (r-1,2),  Q^^\f^)  will  be  approximately 

correct  to  first  order  (as  can  be  seen  from  equation  (23)).  Higher  order  approxi- 
mations, (correct  to  order  p-1  at  f=f^  or  f^),  may  also  be  obtained 

as  in  equation  (22).  It  should  be  noted  that  such  approximations  to  Q(f)  which 
include  the  interaction  between  the  positive  and  negative  frequencies  would  be 
difficult  to  obtain  using  Woodbury’s  identity.  This  is  because  the  identity  would 
have  to  be  applied  to  both  the  positive  and  negative  frequencies  and  would  result 
in  a complicated  expression  for  Q(f)  which  would  provide  little  insight  into  the 
interaction  between  the  various  frequency  components  of  the  sinusoids. 

As  a numerical  example  to  illustrate  the  difference  between  the  approxima- 

1 

tions  for  Q(f)  given  by  (20)  and  (22)  and  the  exact  expression  given  by  (13), 

consider  the  case  when  f = .25;  f„  = .26;  L = 9;  and  For 

12  12  0 

this  case,  the  magnitudes  of  the  coupling  coefficients  between  the  positive  and 
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negative  frequencies  are  (from  (8)):  1 >24  ' " ' >'42  ' " ' ‘>'23  ' " ' "^32  ' " 

I 1 = I I = .033;  1 1 = I | = .034,  and  the  magnitude  of  the  cou- 

pling coefficients  between  the  closely  spaced  lines  isrjy  l=ly  l = ly  1 = 

2.2  22.  34 

I *^43  I = .306.  As  seen  from  equation  (A-4),  equation  (18)  converges  for  this 
case.  In  Figure  1 plots  of  Q(f)  (exact),  Q^°\f),  and  Q^^\f)  are  presented. 

As  is  seen,  Q^^^(f)  provides  a better  approximation  to  Q(f)  than  does  Q^^^(f). 
It  should  be  pointed  out  that  the  value  of  L used  in  this  example  (L  = 9)  was 
only  chosen  for  purposes  of  comparing  the  different  approximations  for  Q(f) 
and  does  not  represent  an  optimal  choice  for  L for  this  example.  Indeed,  as 
previously  discussed,  when  the  correlation  function  is  known  exactly  there  is 
no  cut-off  value  for  L and  the  theoretical  MEM  resolution  improves  without 
limit  as  L is  increased. 

MEM  SPECTRAL  ANALYSIS  OF  MULTIPLE  SINUSOIDS 
IN  1-POLE,  LOW  PASS  NOISE 

In  this  section  we  will  consider  the  case  wiien  the  correlation  lags  are  known 
exactly  and  can  be  written  as 


2 -ajjtl 

p(i)  = CT  e ' ' + 


N 

r 

n=l 


a cos2fff  i 
n n 


(24) 


Equation  (24)  corresponds  to  the  sum  of  N sinusoids  in  1-pole,  low  pass  noise. 
A determination  of  S(f)  for  this  case  indicates  the  resolution  properties  of  MEM 


18 


r 


in  a background  noise  with  variable  spectral  density.  For  this  case  we  expand 
the  prediction  coefficients  as  follows 


2N  2iTjf  k 

g(k)  = C 6(k)  + C„6(k  - L + l)+EAe  ° . (25) 

^ ^ n=l  ° 


where,  as  in  equation  (6),  ^ ~ ^n  particular  choice  for 

the  assumed  solution  is  similar  to  the  solution  of  equation  (4)  for  the  more  gen- 
eral case  when  the  spectral  density  of  the  data  contains  both  poles  and  zeroes 

(Satorius  and  Zeidler,  1977)  and  leads  to  precisely  2N+2  equations  for  the  A , 

n 

C^,  and  Cg*  Substituting  (25)  into  (4)  with  p(je)  given  by  (24)  and  equating 
coefficients  of  exp(2irjf^£),  r=l,2, . . . ,2N,  and  exp(JO!i)  in  the  resulting  equa- 
tion leads  to  the  following  2N+2  equations  for  the  A^,  C^,  and 


2ffjf 


A + 
r 


2N 
S H 
n=l 
n^ 


-C, 


A =- 
rn  n 


-2ffjf^(L-l) 


(26) 


r=l,  2, . . . , 2N 
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wv  2N  a 
„ _ "O!  — n 

" a+2ffjf 

1 - e 


2ffjf  L 
2N  A e ” 

0=7;  S 

^2  , 2Jrjf 


n=l  a 
e -e 


n 


(27) 


(28) 


In  equation  (26)  is  given  by 


^rn  2 - ,2 

L+2a  /a 


(29) 


and  a'^  is  given  by 


2 1 - cos27rf^/cosha  ^ 
tanha 

u 2 _ 2 

where  (n-1, . . . , N).  Substituting  equations  (27)-(28)  into  (26)  leads 

to  the  following  set  of  2N  equations  for  the  A 

r 


2N  / 2Trjf^  \ 

Ap  + S yjTjA^=\e  r=l,2,,..,2N  (31) 

n=l 
n/r 
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where  the  are  given  by 

/3^  =|l  + + (®"^  - cos27rf^)/(cos27rf^  ~ cosha)|  (32) 

and  the  are  given  by 

y'  ) 

* ' *- :SSif — ®r- 

l-e 

Equation  (31)  is  similar  in  structure  to  (7)  and  becomes  identical  to  (7)  as 
O!-*  «,  as  it  should.  Although  the  in  (33)  are  considerably  more  compli- 
cated than  those  for  the  white  noise  case  (equation  (8)),  they  still  approach  zero 
as  L -*  00.  As  the  y^^  approach  zero,  equation  (31)  decouples  and  the  are 
given  to  a good  approximation  by 


As  in  the  case  of  white  noise  (equation  (14)),  Q(f)  evaluated  at  the  frequency  of 

th  2 2 

the  r sinusoid  can  be  simply  expressed  in  terms  of  A and  a /a'1  The  result 

r or 

for  Q(f  ) is 
r 


As  the  approach  zero  and  for  L Q(f^)  is  given  to  a good  approxi- 

mation by  (from  (34)) 

Q(  V 


r 


Equation  (38)  shows  that  for  large  L,  the  values  Q(f^)  are  a funct  ion  of  f 
2 2 

as  well  as  o /a  . This  result  is  reasonable  since  the  signal -to-noise  ratio 
r o 

(SNR)  per  unit  bandwidth  is  a function  of  frequency.  This  indicates  that  the 
resolution  capabilities  of  the  MEM  spectral  estimate  will  be  dependent  on  where 
the  sinusoids  are  located  in  the  noise  spectrum.  For  sinusoids  located  in  the 
low  SNR  regions,  the  MEM  resolution  will  be  worse  than  for  sinusoids  located 
in  the  high  SNR  region.  This  is  indicated  in  Figure  2 where  Q(f)  is  plotted  for 
4 sinusoids  of  equal  amplitudes.  Two  of  the  sinusoids  are  located  near  the  SNR 
region  and  the  other  two  sinusoids  are  located  near  the  high  SNR  region.  As  is 
seen,  the  two  sinusoids  in  the  low  SNR  region  are  poorly  resolved  whereas  the 
two  sinusoids  in  the  high  SNR  region  are  well  resolved  (arrows  in  Figure  2 in- 
dicate the  correct  location  of  the  frequencies  of  the  sinusoids).  As  in  the  case 
of  Figure  1,  the  value  of  L used  in  Figure  2 (L=14)  was  only  chosen  for  pur- 
poses  of  illustrating  the  variation  of  two  sinusoid  MEM  resolution  versus  fre-  i 

quency  and  does  not  represent  an  optimal  value  for  L. 

It  should  be  noted  that  in  complete  analogy  with  the  white  noise  case,  the 
proper  time  series  model  for  N sinusoids  in  additive  1-pole,  low  pass  noise  is 
an  ARMA  model  which  contains  2N+1  AR  terms  and  2N  MA  terms.  (The  ' 

additional  AR  term  is  due  to  the  1-pole  structure  of  the  additive  noise. ) There- 
fore, one  expects  that  as  L ■*  ®,  S(f)  will  converge  to  the  power  spectral  den- 
sity of  the  1-pole  noise  (for  f f ; n=l, . . . ,N).  This  can  easily  be  seen  from 

n 

the  development  presented  in  this  section  by  noting  from  equation  (34)  that 
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Ito  (f  -f)  - 0. 
L”*“ 


provided  f ^ f (n=l» . . • ,N).  Further  from  equations  (27)  and  (28),  it  is  seen 
_2 

that  C “*  c and  C "•  0 as  L "•  ®.  Therefore,  from  equations  (2),  (35),  and 

X ^ 

(36)  we  have 


lim  S(f) 

L-oo 


A 2 

^ sinha 
o 

coshtt  - cos2Pf  ’ 


f ^ f^(n=l,...,N)  . 


(39) 


Equation  (39)  is  the  expression  corresponding  to  the  power  spectral  density  of  the 
1 -pole,  low  pass  noise.  Therefore,  from  equations  (38)  and  (39)  it  is  seen  that 
S(f)  provides  an  increasingly  accurate  pointwise  approximation  to  the  background 
noise  spectrum  as  well  as  precise  sinusoid  resolution  capabilities  as  L -•  ®. 
However,  as  noted  in  the  previous  section,  this  only  applies  when  the  correla- 
tion lags  are  known  exactly. 
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CONCLUSIONS 


In  this  papor,  the  theoretical  MEM  spectral  estimate  of  multiple  sinusoids 
in  noise  has  been  examined  by  the  method  of  undetermined  coefficients.  For 
the  case  of  sinusoids  in  white  noise,  the  prediction  filter  coefficients  were  ex- 
panded directly  in  terms  of  the  input  sinusoids  (equation  (6)).  This  expansion 
leads  to  a transformation  of  the  original  L X L equations  for  the  prediction 
coefficients  (equation  (4))  to  a set  of  2N  x 2N  equations  (equation  (7)).  The 
transformed  equations  are  a smaller  set  of  linear  coupled  equations  than  equa- 
tion (4)  when  L > 2N  and  are  particularly  useful  for  purposes  of  computing  the 
theoretical  MEM  spectral  estimate  for  large  L (high  resolution  limit).  Further, 
the  reduced  equations  (equation  (7))  have  been  shown  to  provide  additional  insight 
into  the  analytical  structure  of  the  MEM  spectral  estimate.  Also,  for  certain 
cases  where  there  is  little  interaction  between  some  of  the  frequency  compc- 
nents  of  the  input  sinusoids,  the  2N  x 2N  equations  may  be  approximated  by  a 
smaller  set  of  coupled  equations. 

For  the  case  of  N sinusoids  in  1-pole,  low  pass  noise,  the  prediction  filter 
coefficients  were  expanded  in  terms  of  the  input  sinusoids  as  well  as  two  delta 
functions  which  occur  at  the  beginning  and  end  of  the  prediction  filter  (equa- 
tion (25)).  This  expansion  also  leads  to  a set  of  2N  x 2N  equations  (equation 
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(31)).  For  this  case,  the  values  of  the  MEM  spectral  estimate  evaluated  at  the 
input  sinusoid  frequencies  were  shown  to  be  a function  of  the  frequency  of  the 
sinusoids.  This  result  is  reasonable  since  the  SNR  per  unit  bandwidth  is  also 
a function  of  the  frequency. 

The  results  derived  in  this  paper  give  a further  understanding  of  the  MEM 
spectral  estimate  of  sinusoids  in  noise  when  the  autocorrelation  lags  are  esti- 
mated from  the  data.  As  the  variance  in  the  estimates  of  the  lags  becomes 
appreciable,  deviations  in  the  MEM  spectral  estimate  wUl  being  to  appear, 
especially  at  its  peaks  (e.g. , Lacoss,  1971,  Baggeroer,  1976). 
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Figure  2. 


A plot  of  Q(f)  for  4 sinusoids  in  I pole  low  pass  noise. 
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APPENDIX  A.  CONDITIONS  FOR  CONVERGENCE  OF  PERTURBATION 
EXPANSION  IN  EQUATION  (18) 

It  is  the  purpose  of  this  Appendix  to  establish  sufficient  conditions  for  the 
convergence  of  the  perturbation  expansion  expressed  by  equation  (18).  We  make 
use  of  the  following  basic  result  (e.g. , Stewart,  1973):  The  perturbation  series 
I + P + . . converges  to  Q - P)  ^ if  1 Ig  1 1 < 1»  where  11*11  denotes  any 

matrix  norm  such  that  1 1a  *C  1 1 S i 1 A 1 1 1 1 C 1 1 . 

The  particular  norm  of  which  we  make  use  and  which  is  relatively  simple  to 
calculate  is  the  “ - norm  which  is  denoted  by  11*11^.  The  » - norm  of  a 
p X p matrix,  A,  is  defined  as  follows  (Stewart,  1973) 


where  A.,  are  the  elements  of  A.  Therefore,  a sufficient  condition  for  the 
1] 

convergence  of  equation  (18)  is 


Ti 
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Since  1 1 AC  ^ I A 1 Ic  | | , then  a simpler  sufficiait  condition  to  check  is 

ssszOO  sOD^OB 


(A-3) 


For  the  special  case  of  two  sinusoids  of  equal  power  in  white  noise  in  which 
there  is  little  interaction  between  the  positive  and  negative  frequency  components 
of  the  sinusoids,  condition  (A-3)  reduces  to 


1 - ly, 


max  + Iv  J,  |y„,l  + ly„,  1 


12 


'13 


'14 


23' 


'24 


<1 


(A-4) 


where  the 


are  computed  from  equation  (8). 
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